Explosive first order transition to synchrony in networked chaotic oscillators 



o 

(N 



(N 



I. Leyva,^'^ R. Sevilla-Escoboza,^ J. M. Buldu/'^ I. Sendifia-Nadal,^'^ J. Gomez-Gardefies,^' ^ 

A. Arenas,^' ^ Y. Moreno,^' ^ S. Gomez, ^ R. Jaimes-Reategui,^ and S. Boccaletti^ 

^Complex Systems Group, Univ. Rey Juan Carlos, 28933 Mostoles, Madrid, Spain 

^Center for Biomedical Technology, Univ. Politecnica de Madrid, 28223 Pozuelo de Alarcon, Madrid, Spain 

^ Centro Universitario de Los Lagos, Univ. de Guadalajara, Lagos de Moreno, Jalisco 4^460, Mexico 

^Dep. de Fisica de la Materia Gondensada, Univ. de Zaragoza, Zaragoza 50009, Spain 

^Institute for Biocomputation and Physics of Gomplex Systems (BIFI), Univ. de Zaragoza, 50009 Zaragoza, Spain 

^Dep. d'Enginyeria Informdtica i Matemdtiques, Univ. Rovira i Virgili, 43007 Tarragona, Spain 

^Dep. de Fisica Teorica, Univ. de Zaragoza, Zaragoza 50009, Spain 

Critical phenomena in complex networks, and the emergence of dynamical abrupt transitions in 
the macroscopic state of the system are currently a subject of the outmost interest. We report 
evidence of an explosive phase synchronization in networks of chaotic units. Namely, by means of 
both extensive simulations of networks made up of chaotic units, and validation with an experiment 
of electronic circuits in a star configuration, we demonstrate the existence of a first order transition 
towards synchronization of the phases of the networked units. Our findings constitute the first prove 
of this kind of synchronization in practice, thus opening the path to its use in real- world applications. 
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The understanding of the spontaneous emergence of 
collective behavior in ensembles of networked dynamical 
units constitutes a fascinating challenge in science. De- 
spite the fact that critical phenomena in networks have 
been intensively studied, the physics literature [1 almost 
exclusively reports continuous phase transitions. How- 
ever, it has been recognized that, although very few in 
number [2tt4j, there are physical processes which might 
lead to sharp, discontinuous transitions of a global order 
parameter. The last several years have also witnessed 
an ever-increasing interest in studying networked sys- 
tems composed of nonlinear dynamical units f5], and in 
particular, in the emergence of synchronization phenom- 
ena [6 . Within this latter context, some advances have 
been made for the case of non-equilibrium synchroniza- 
tion transitions of chaotic systems [7l[8], being, however, 
all the reported cases examples of second order phase 
transitions. 

More recently it has been shown that discontinuous 
transitions can take place in networks of periodic os- 
cillators [4 . There, a first-order non-equilibrium syn- 
chronization transition has been proved in scale-free net- 
works, as a consequence of a positive correlation between 
the heterogeneity of the connections and the natural fre- 
quencies of the oscillators. Whether such an explosive be- 
havior is restricted to periodic oscillators, or can be gen- 
eralized to more complicated dynamical units remained 
an open problem. 

In this Letter, we investigate the critical properties 
of synchronization transitions in heterogeneous networks 
of chaotic units (when the aforementioned interplay be- 
tween topology and dynamics is taken into account), and 
we will give the first numerical and experimental evidence 
of an explosive phase synchronization in such ensembles 
of complex systems. Even though the extent of our dis- 
cussion is valid regardless of the specific phase coherent 



[9^ dynamics that is considered for the evolution of each 
network's node, from now on we will focus on a specific 
benchmark chaotic system. Additionally, our choice is 
dictated by the need of implementing a qualitatively sim- 
ilar model in the laboratory, and therefore by the unique 
opportunity of contrasting the numerical predictions with 
the experimental evidence. 

Let us consider an ensemble of A^ = 1,000 piecewise 
Rossler units, interacting in a network via a bidirectional 
diffusive- like coupling [TOl [11] : 
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Every node (indexed by i = 1,...,A/') is here repre- 
sented by an associated three-dimensional vector Xi(t) = 
{xi(t)^yi(t)^Zi(t)). The parameters are: F = 0.05, j3 = 
0.5, A = 1, /i = 15, and v = 0.02 — ^, where R is a tun- 
able quantity that regulates the dynamical state of the 
system. In particular, R induces a chaotic dynamics [J^ 
in the range i? = [55, 110]. A = {a^j} is the adjacency 
matrix (a^j = 1 if units i and j are connected, and 
otherwise), and d is the coupling strength. Notice that 
the coupling is here applied through the x variable, but 
equivalent results are obtained using a coupling in the 
variable y. 

Finally, the natural oscillation frequency of node i de- 
pends linearly on the parameter a^. Therefore, as long 
as one is concerned with imposing a positive correlation 
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FIG. 1: (Color online). Phase synchronization degree 5* as a 
function of the coupling strength d for different SF networks of 
size A^ = 1, 000, and average degree {k) = 6. The networks are 
built using the configuration model [14] (see text for details). 
7 = 2.2 (red triangles), 7 = 2.5 (green squares), 7 = 3.0 (blue 
diamonds and black circles). The correlation between node 
degree and natural frequency is set via Eq. ([3|. Aa = 6.0 
and R — 70 for the black circles case, while Aa = 10.0 and 
i? = 100 for the other networks. Continuous (dashed) lines 
mark the forward (backward) simulations, as d is increased 
(decreased) in steps of dd = 3 x 10~^. The inset plot reports 
the average synchronization error (e) vs. d in the proximity of 
the first-order transition occurring for 7 = 3.0, R = 100 and 
Aa = 8.0. 



between the natural frequency and the degree ki of each 
unit (the number of connections that the i^^ unit is form- 
ing with the rest of the network), the a^ values (and 
therefore, the oscillators' frequencies) are distributed fol- 
lowing the relation: 
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where a = 10^, and Aa is a factor that determines the 
slope of the linear distribution. As a consequence, all 
nodes with degree k = 1 have the same natural frequency 
regardless of the value of Aa. On the other hand, the 
range of frequencies in the ensemble becomes wider as 
Aa and /cmax (the maximum degree in the network) are 
increased, i.e. as more heterogeneous degree distributions 
and/or steeper slopes are considered. 

The state of the network is monitored as a function 
of the coupling, by gradually increasing the value of d 
in steps dd along the simulation from d = (without 
resetting the system, as it will be done later in the ex- 
periment). Whenever the coupling is increased in Sd, a 
long transient is discarded before the data are used for 
further processing. Furthermore, as long as we are look- 
ing for a first-order phase transition (and thus for an 
expected associated synchronization hysteresis), we per- 
form the simulations also in the reverse way, i.e. starting 
from a given value (imax where the ensemble is phase syn- 
chronized, and gradually decreasing the coupling by Sd 



at each step. In what follows, the two sets of numerical 
trials are termed as forward and backward^ respectively. 
We here focus in the phase synchronization regime 
[9l [12]. For each node z, the instantaneous phase 
at time t is geometrically evaluated [T3] as ^^(t) = 
arctan {yi{t)/xi{t)), and the mean synchronization degree 



is calculated as S 
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tical bars denote the module and the angle brackets a 
temporal averaging. 

The first result is that the delicate equilibrium between 
the network heterogeneity and the frequency distribution 
plays a crucial role in determining the nature of the phase 
synchronization transition. This is shown in Fig. 1, where 
we compare the dependence of the synchronization degree 
S on the coupling strength d, in different networks of 
same size and average degree. All graphs are scale- free 
(SF) networks, obtained with the configuration model 
^14] , and featuring a degree distribution P{k) ~ A:~^, 
with 7 = 2.2 (red triangles), 7 = 2.5 (green squares), and 
7 = 3.0 (blue diamonds and black circles). Remarkably, 
the specific chaotic state of the nodes is a relevant feature, 
as demonstrated by comparing the case Aa = 6.0, i? = 70 
(black circles, second order phase transition) with the 
case Aa = 10.0, i? = 100 (blue diamonds, first order 
phase transition) for 7 = 3.0. The explosive (first order) 
character of the transition is clearly manifested by the 
presence of a hysteresis region, delimited in the figure by 
the continuous and dashed lines that mark, respectively, 
the forward and backward cycles. The inset plot shows 
the average synchronization error (e) (the average of the 
Euclidean distance between all pairs of nodes, normalized 
to account for the dependence of the mean amplitude 
oscillations on the parameters) vs. d, in the proximity of 
the transition for 7 = 3.0, R = 100 and Aa = 8.0. 

The conclusive information conveyed by Fig. 1 is that 
the collective behavior emerging in the system is, in- 
deed, a synchronization of the phases of the chaotic units, 
which is not associated to any specific amplitude correla- 
tion (as the error (e) is even increasing at the right side 
of the transition point). 

A more exhaustive description of the phenomenon 
emerges from the exploration of the nature of the tran- 
sition in the full parameter space d — R. The results are 
shown in Fig. 2, where we report (in color code) the val- 
ues of S^ for both forward and backward simulations. As 
it can be seen, the plane d — R can be clearly divided in 
two areas (denoted as / and // in the figure) where the 
transition is of the first and second order, respectively. 
The striped portion of the area where the transition is 
of the first order marks the region where the hysteresis 
phenomenon is observed. 

A first conclusion of this Letter is that our numeri- 
cal study reveals the presence, as a generic feature, of 
a first-order non-equilibrium phase transition towards 
phase synchronization in a network of chaotic oscilla- 



100 


m^.. 


90 


Wmr 


80 


fc 


70 


^^^k 



0.5 



9 10 



d(xia^) 



FIG. 2: (Color online). Mean synchronization degree S (see 
text for definition) in the parameter space d — R (the values 
of S are color-encoded, according to color-bar at the right of 
the panel). The horizontal dashed line marks the separation 
between the region of the parameter space where a second- 
order transition occurs (below the line), and that for which 
the transition is instead of the first-order type (above the 
line). In this latter region, a typical hysteresis phenomenon 
in the forward and backward simulations is observed (whose 
extent is here marked by a striped area). The step used for 
the variations in the coupling strength is (5c? = 5 x 10~^. The 
entire phase diagram refers to the case of a scale free network 
with A^ = 1, 000, (/c) = 6, 7 = 3 and Aa = 6.0. 



tors. The phenomenon is here originated by a setting 
of the constituent nodes within a proper chaotic regime, 
together with imposing a positive correlation between the 
heterogeneity of the graph connections and the natural 
frequencies of the oscillators. It is worth noticing that in 
principle one can have first-order transitions for more ho- 
mogeneous topologies, however at the cost of expanding 
the range of natural frequencies in the network. 

We now move to validate the robustness and general- 
ity of our results against stochastic fluctuations, and/or 
parameter mismatches. For that purpose, we implement 
the electronic network schematically shown in Fig. 3. The 
experiment consists of six piecewise Rossler circuits oper- 
ating in the chaotic regime, which are labeled as A/"!, A/'2, 
. . . , NQ. The details of the circuit construction, as well 
as the qualitative equivalence of the experiment with the 
model of Eq. (1) are available in the literature [TOlllI]. 

The circuits are arranged in a star-like configuration, 
which, on its turn, represents the maximally heteroge- 
neous structure available for small ensembles. All the 
chaotic oscillators have the same internal i^exp [the ex- 
perimental equivalent of the parameter R in Eq. (1)], 
to ensure that they work in an almost equal dynamical 
regime. By a fine tuning of the values of the capacitors, 
circuits are configured such that the central node Nl os- 
cillates with a mean frequency of 3, 333 Hz, and the leaves 
nodes A/'2, . . . , NQ are set with frequencies in the range 
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FIG. 3: (Color online). Sketch of the experimental setup. Six 
Rossler circuits (blue nodes) are connected in a star config- 
uration. The bidirectional coupling is adjusted by means of 
five digital potentiometers X9C104 (red nodes) whose param- 
eters Cu/d (value of the resistance) and Cstep (increment of 
the resistance at each step) are controlled by a digital signal 
coming from a DAQ Card. The outputs of the six circuits 
are sent to a set of voltage followers that act as a buffer and, 
then, sent to the analog ports (AI 0, AI 1, . . . , AI 5) of the 
same DAQ Card. The whole experiment is controlled from a 
PC with a Labview Software. 



of 2, 240 ± 200 Hz. Notice that, due to the experimental 
variability, the frequencies of the oscillators suffer also 
from unavoidable dispersions. 

The experiment is performed by coupling simultane- 
ously the six oscillators (blue nodes) with the same cou- 
pling parameter d. We use five digital potentiometers 
(red nodes) in order to warrant that the parameter d 
is changed simultaneously for all nodes. Digital poten- 
tiometers (X9C104) are adjusted by a digital signal com- 
ing from ports PO.O and PO.l of a NI Instruments DAQ 
Card (NI DAQPad-6259). The output of each circuit is 
connected to a voltage follower that works as a buffer. 
Next, signal is acquired by the analog ports (AI 0, AI 1, 
. . . , AI 5) of the same DAQ Card, and recorded in a 
PC for further analysis. The incoming signal of the ana- 
log inputs (ADC) and the signal sent through the digital 
outputs (DOP) are controlled and recorded by a Labview 
Software. 

Once the data are stored, processing of them is made 
to obtain the equivalent instantaneous phases. These 
latter quantities are calculated by defining the instan- 
taneous phase 0i(t) of each oscillator i as (j)i{t) = 27r/^ + 



in each interval ti. < t < t^.+i, where ti. 



IS 



tii+i—ti^ 

the instant at which the Ij^ crossing of the i^^ oscillator 
with its Poincare section occurs. It is crucial to remark 
that, for phase coherent systems, such a measure of in- 
stantaneous phase is fully equivalent to the geometrical 
phases used in the numerical data [9 . 

In Fig. 4 we report the experimental synchronization 
diagram for three values of the parameter i^exp, both for 
forward and backward variations of the coupling strength 
d. As in the numerical case, the dynamical regime tuned 
by Rexp determines whether the transition is of the first 
or second order. The former case can be identified by the 
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FIG. 4: (Color online). Phase synchronization degree S vs. d 
for the star configuration shown in Fig. 3. The values of the in- 
ternal resistance determining the chaotic state are Rexp = 55 
(black circles, second-order phase transition), i?exp = 65 (red 
squares, first order phase transition with narrow hysteresis) 
and Rexp = 70 (blue stars, first order phase transition with 
wide hysteresis). The hysteresis in the forward (red and blue 
solid lines) and backward (red and blue dashed lines) exper- 
imental trials reveals the genuine character of a first-order 
phase transition. The coupling strength varies at steps of 
7.5 X lO"'^ and the values of the oscillator frequencies are 
given in the text. 



not need to fine-tune the coupling strength to get the 
transition, but to vary the system parameters {R in our 
case) in a wider region. These two observations are of 
utmost importance when it comes to translate the un- 
covered mechanism into practice. We then expect that 
our work will open the path to the use of explosive syn- 
chronization phenomena in many relevant applications. 
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sharp increase of the synchronization parameter 5, and 
the existence of a hysteresis. Interestingly, as can be seen 
in Fig. 2 for the two larger values of i^exp, this parameter 
also determines the range of d at which the hysteresis 
appears (i.e. the width of the hysteresis window). Finally, 
it has to be noticed that the values of i?exp at which 
the order of the transition changes display mismatches 
with those of the numerical trials, due to the imperfect 
equivalence between model and experimental setup. 

In conclusion, we have given here numerical and exper- 
imental proof of the emergence of a first order synchro- 
nization transition in a network of phase coherent chaotic 
oscillators. Such a discovery is the first evidence of this 
new paradigm of synchronization in chaotic systems. We 
have shown that the correlation between oscillation fre- 
quencies and structural heterogeneity is a necessary con- 
dition for such an explosive phenomenon. Furthermore, 
our work goes one step further by also establishing that, 
in chaotic systems, additional features like the specific 
dynamical state at which the chaotic units operate are 
determinant. This is of particular interest mainly for 
two reasons: first, real systems are noisy, and do not al- 
ways work on the same regime. Therefore, whether or not 
a first order synchronization transition can be achieved 
would not depend on the specifics of the system, but on 
the region at which it is operating. Secondly, one does 



[1 

[2: 
[3: 

[5: 
[6: 

[r 

[9: 
[lo: 

[11 

[12: 
[13: 

[14 



S. N. Dorogovtsev, A. V. Goltsev and J. F. F. Mendes. 

Rev. Mod. Phys. 80, 1275-1335 (2008). 

P. Echenique, J. Gomez- Gardehes and Y. Moreno. Euro- 

phys. Lett. 71, 325 (2005). 

M. Karsai, J-Ch. Angles d'Auriac, and F. Igloi. Phys. 

Rev. E 76, 041107 (2007). 

J. Gomez- Gardehes, S. Gomez, A. Arenas and Y. 

Moreno. Phys. Rev. Lett. 106, 128701 (2011). 

S. Boccaletti, V. Latora, Y. Moreno, M. Chavez and D. 

U. Hwang. Phys. Rep. 424, 175-308 (2006). 

A. Arenas, A. Diaz-Guilera, J. Kurths, Y. Moreno and 

C. S. Zhou. Phys. Rep. 469, 93-153 (2008). 

C. J. Tessone, M. Cencini and A. Torcini. Phys. Rev. 

Lett. 22, 224101 (2006). 

M. Cencini, C. J. Tessone and A. Torcini. Chaos 18(3), 

037125 (2008). 

S. Boccaletti, J. Kurths, G. Osipov, D. L. Valladares and 

C. S. Zhou. Phys. Rep. 366, 1-101 (2002). 

A. N. Pisarchik, R. Jaimes-Reategui, J. R. Villalobos- 

Salazar, J. H. Garcia-Lopez and S. Boccaletti. Phys. Rev. 

Lett. 96, 244102 (2006). 

A. N. Pisarchik, R. Jaimes-Reategui, R. Sevilla-Escoboza 

and S. Boccaletti. Phys. Rev. E 79, 055202(R) (2009). 

M. G. Rosenblum, A. S. Pikovsky and J. Kurths. Phys. 

Rev. Lett. 76, 1804-1807 (1996). 

G. V. Osipov, A. S. Pikovsky, M. G. Rosenblum and J. 

Kurths. Phys. Rev. E 55, 2353-2361 (1997). 

E. A. Bender and E. R. Canfield. J. Comb. Theory A 24, 

296-307 (1978). 



